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Currently the best prospect for detecting Weakly Interacting Massive Particles (WIMPs) is via 
the annual modulation, which occurs due to the Earth's rotation around the Sun, of the direct 
detection signal. We investigate the effect of uncertainties in our knowledge of the structure of the 
galactic halo on the WIMP annual modulation signal. We compare the signal for three non-standard 
halo models: Evans' power-law halos, Michie models with an asymmetric velocity distribution and 
Maxwellian halos with bulk rotation. We then compare the theoretical predictions of these models 
with the experimental signal found by the DAMA experiment and investigate how the WIMP mass 
and interaction cross section determined depend on the halo model assumed. We find that the 
WIMP mass confidence limits are significantly extended to larger masses, with the shape of the 
allowed region in the mass-cross section plane depending on the model. 
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I. INTRODUCTION 

The rotation curves of spiral galaxies are typically flat 
out to about ~ 30 kpc. This implies that the mass en- 
closed increases linearly with radius, with a halo of dark 
matter extending beyond the luminous matter (!]]. The 
nature of the dark matter is unknown with possi- 
ble candidates including MAssive Compact Halo Objects 
(MACHOs), such as brown dwarves, Jupiters or black 
holes and elementary particles, known as Weakly Inter- 
acting Massive Particles (WIMPs), such as axions and 
the neutralino, the lightest supersymmetric particle. 

WIMPs can be directly detected via their elastic scat- 
tering off target nuclei. In the long term the directional 
dependence of the detector recoil will provide the best 
means of direct WIMP detection ||. Currently the best 
prospect for distinguishing a WIMP signal from the de- 
tector background is via the annual modulation of the 
signal, which occurs due to the Earth's rotation around 
the sun The DAMA collaboration f|, using a de- 

tector consisting of radiopure Nal crystal scintillators at 
the Gran Sasso Laboratory, have recently reported a 4tr 
annual modulation signal consistent with the detection 
of a WIMP with mass m v = 52±|° GeV @, assuming a 
Maxwellian halo with velocity dispersion roughly equal 
to the local rotation velocity, vq = 220kms~ 1 . 

The values of the WIMP mass and interaction cross 
section found from the annual modulation signal are 
known to depend somewhat on the assumed values of 
poorly known astrophysical parameters. The effects 
of bulk rotation [^|j9| and varying the velocity disper- 
sion jlO|,[)) , within the observationally allowed range vo = 
220±40kms _1 |flj], have been examined for a Maxwellian 
halo. The DAMA collaboration have subsequently in- 
cluded these uncertainties in the analysis of their latest 
data 0. 

The standard Maxwellian halo model has a number 
of deficiencies (for example see Ref. [n2[ and references 



therein). The halo may not be spherical; N body sim- 
ulations of gravitational collapse produce axisymmetric 
or triaxial halos, and indeed other spiral galaxies appear 
to have flattened halos (p|. The power-law halo mod- 
els of Evans |l4|] provide an analytically tractable frame- 
work for investigating the effect of varying halo prop- 
erties such as the flattening, and have previously been 
used to investigate variations in the mean total |15| and 
directional [|| WIMP detection rates. Another possible 
modification to the standard halo model is an asymmet- 
ric velocity distribution and the annual modulation 
of the total @@Jl|] and directional fu| WIMP signals 
have been calculated for various asymmetric velocity dis- 
tributions. 

In this paper we compare the variation of the annual 
modulation signal for power-law halo models, asymmet- 
ric velocity distributions and Maxwellian halos with bulk 
rotation. We then examine how the range of WIMP 
masses consistent with the latest DAMA data depends 
on the halo model assumed. 



II. ANNUAL MODULATION SIGNAL 

The WIMP detection rate depends on the speed dis- 
tribution of the WIMPs in the rest frame of the detector, 
/„. This is found from the halo velocity distribution, /(v) 
by making a Galilean transformation v — * v = v + v , 
where v is the Earth's velocity relative to the galactic 
rest frame, and then integrating over the angular distri- 
bution In galactic co-ordinates the axis of the ecliptic 
lies very close to the <p — z plane and is inclined at an 
angle 7 ~ 29.80° to the <fi — r plane Eq|. Including all 



"The effect of the Sun's gravity on the WIMP distribution 
can be neglected |l8|. 
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components of the Earth's motion, not just that parallel 
to the galactic rotation: 

v c = V\ sin af + 

(vq + v\ cos a sin 7) <f> — v\ cos a cos 7 z , (1) 

where vo ~ 232kms~ 1 is the speed of the sun with respect 
to the galactic rest frame, v\ w 30kms _1 is the orbital 
speed of the Earth around the Sun and a = 2Tr(t~to) /T, 
with T — 1 year and to ~ 153 days (June 2nd), when the 
component of the Earth's velocity parallel to the Sun's 
motion is largest. 

In the range of masses and interaction cross sections ac- 
cessible to current direct detection experiments the best 
motivated WIMP candidate is the neutralino, for which 
the event rate is dominated by the scalar contribution. 
The differential event rate can then be written in terms 
of the WIMP cross section on the proton, 
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where m p is the proton mass and f p the effective WIMP 
cross section on the proton. The differential event rate 
simplifies to (see the Appendix for more details): 



dR ^ 
dE=^ 
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where the local WIMP density, p x has been normalised 
to a fiducial value po.3 = 0.3GeV/cm 3 , such that £ = 
Px/Po.3, E is the recoil energy of the detector nucleus 
and T(E) is defined as p| 
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— dv . 
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where 



is the minimum detectable WIMP velocity 
'E(m x +m A )^ 1/2 



2m x mA 



(5) 



and rriA the atomic mass of the target nuclei. 

In order to compare the theoretical signal with that 
observed we need to take into account the response of the 
detector. The electron equivalent energy, E oc , which is 
actually measured is a fixed fraction of the recoil energy: 
-Eoc = QaE. The quenching factors for I and Na are 
qi = 0.09 and </Na = 0.30 respectively pof . The energy 
resolution of the detector || is already taken into account 
in the data released by the DAMA collaboration. 

The expected experimental spectrum per energy bin 
for the DAMA collaboration set-up is then given by |1G] 
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Energy (keV) 


S ,k (cpd/kg/keV) 


S m , k (cpd/kg/keV) 


2-3 


0.54 ± 0.09 


0.023 ± 0.006 


3-4 


0.21 ± 0.05 


0.013 ± 0.002 


4-5 


0.08 ± 0.02 


0.007 ± 0.001 


5-6 


0.03 ± 0.01 


0.003 ± 0.001 



TABLE I. So,k and S m , k values obtained, by the DAMA 
collaboration, from a maximum likelihood analysis of their 
combined data from four annual cycles. 



where rN a = 0.153 and r\ — 0.847 are the mass fractions 
of Na and I respectively. 

Since vo S> v\ the differential event rate in the k-th en- 
ergy bin can be expanded in a Taylor series in cos a [ElJ : 



AR 
~AE 
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(7) 



The DAMA collaboration use a maximum likelihood 
method to separate the time-independent background 
from the WIMP signal and have released the resulting 
values of So k and S'm.k, and the errors on them, for each 
energy bin M (see Table 1). 



III. HALO MODELS 

Since the resolution of numerical simulations is not yet 
large enough to allow the numerical calculation of the an- 
nual modulation signal |22|, we have to use analytic mod- 
els for the velocity distribution of the dark matter halo. If 
the galactic halo contains large amounts of substructure, 
as N-body simulations appear to indicate then the 
WIMP detection rate may be dramatically altered if we 
are currently passing through a clump of substructure. 
The extent to which substructure is actually present in 
the halo is currently a matter of debate however. 

The effect of the halo model on the annual modulation 
signal can be assessed most simply, with least recourse 
to the detector properties, in terms of the dimensionlcss 
function T(E) f, as defined in Eq (|J). The yearly aver- 
aged value oiT(E) depends weakly on the halo model |TJ] 
decreasing with increasing recoil energy, with the de- 
crease being most rapid for small WIMP masses. For 
each model we plot the annual variation in T(E): 



AT{E) = 
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(8) 



defined such that AT(E) is taken to be positive if T(E) 
is largest in June (a = 0) as found by the DAMA collab- 
oration, as a function of recoil energy for a Ge 76 detector 



^The shape of the energy spectrum measured also depends 
on the nuclear form factor F(q) however. 
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FIG. 1. The annual variation in T(E), AT(E), for 
halo models with asymmetric velocity distributions with 
A = 0,0.5, 1, from top to bottom row, for four values of the 
WIMP mass m x = 30, 50, 100, 200 GeV, from bottom to top 
for a Ge 76 detector. 



and four values of the WIMP mass: 



= 30,50,100,200 



GeV. Values for other monatomic detectors can be found 
by rescaling the x-axis by tua/ (tua + m x) 2 ■ 



A. Asymmetric velocity distribution 

An extension of the Michie model can be used as a 
simple model of velocity anisotropy p|,|l6|] 



/(v) = V 
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where v csc is the escape velocity and V is a normalisation 
constant. The deviation of the velocity distribution from 
isotropic is parameterised by A, the standard Maxwellian 
halo, cut-off at the escape velocity, is recovered for A = 0. 
Since the halo is formed by gravitational collapse the 
deviation is most likely to be towards radial orbits, with 
< A < 1 |. We take v csc = 650kms _1 , although the 
effect of variations in the value of the escape velocity is 
small §§. 

The annual variation in T(E), AT(E), is plotted in 
Fig. g for A = 0, 0.5, 1. The mean value of T{E) is slightly 
higher for the asymmetric models than for the standard 
Maxwellian halo model. 



B. Power- law halos 

Evans' family of axisymmetric distribution func- 
tions pi lead to velocity distribution functions, in the 
rest frame of the galaxy, of the form JTa] : 
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FIG. 2. The annual variation in T(E), AT(E), for 
power-law halos with q = 1,0.85,0.707, from top to bottom 
row, for four values of the WIMP mass m x = 30, 50, 100, 200 
GeV, from bottom to top for a Ge 76 detector. 
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with 
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where R c is the core radius, Rq is the solar radius and 
q is a flattening parameter, which varies between 1, for 
spherical halo, and 1/V2 = 0.707. Following Ref. |15j we 
take R c = 8.5kpc and Rq — 7kpc and explore the effect 
of varying q. 

The annual variation in T(E), AT(E), is plotted in 
Fig. | for q = 1,0.85,0.707. The mean value of T{E) 
is slightly lower for the flattened halos than for the 
Maxwellian halo. The change in the annual variation 
due to flattening is larger than that due to an asymmet- 
ric velocity distribution, particularly so for large recoil 
energies and small WIMP masses. 



C. Bulk rotation 

Halo models with bulk rotation can be constructed by 
taking linear combinations of the velocity distribution 
function Ittal: 



/rot(v) = a rot /+(v) + (1 - a rot )/_(v) 



(12) 



where 
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FIG. 3. The annual variation in T(E), AT(E), for ro- 
tating Maxwellian halo models with a rot = 0.36, 0.5, 0.64, 
from top to bottom row, for four values of the WIMP mass 
m x = 30, 50, 100, 200 GeV, from bottom to top for a Ge 76 
detector. 
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and a ro t is related to the dimcnsionless galactic angu- 
lar momentum, A rot : A rot = 0.36|a rot — 0.5|. Numerical 
studies of galaxy formation find that |A ro t| < 0.05 
corresponding to 0.36 < a rot < 0.64. A non-rotating halo 
has a ro t = 0.5, whilst a counter-rotating (co-rotating) has 
a rot < 0.5 (a ro t > 0.5). 

The annual variation in T(E), AT(E), is plotted in 
Fig. |3| for Maxwellian halos with a rot = 0.36,0.5,0.64. 
The mean value for counter(co)-rotation is lower (higher) 
than for the non-rotating Maxwellian halo at small recoil 
energies, and higher (lower) at large recoil energies. The 
change in the annual variation, which is suppressed (en- 
larged) for counter(co)-rotation, is largest for small recoil 
energies. 

For each of the halo models studied the annual varia- 
tion is negative for small recoil energies (E < 4 keV for 
M x = 30GeV, E < 30 keV for M x = 200 GeV) increas- 
ing as the recoil energy increases. For fixed recoil energy 
the annual variation is largest for small WIMP masses. 



IV. ANALYSIS OF THE DAMA DATA 

Brhlik and Roszkowski [jl0| have devised a tech- 
nique for comparing the experimental results released by 
DAMA with theoretical predictions for the annual modu- 
lation signal, in the absence of detailed information about 
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FIG. 4. The n = 35 contour, delineating the region of 
m x — £<t parameter space compatible with the DAMA an- 
nual modulation signal, for asymmetric halo models with 
A = 0, 0.5, 1 (solid, dashed and dotted lines respectively). 



the experimental set-up, such as the efficiency of each Nal 
crystal. They define a function n: 



tb 



_ c cx P\2 



(15) 



y 0,k 



where the experimental errors on the time dependent and 
independent parts of the signal, cr m ,k and oo,k respec- 
tively, serve as weights. The contour, in the m x — £er p 
plane, k = 35 agrees reasonably well with the DAMA col- 
laborations 3<7 contour. Whilst this approach does not 
give accurate confidence limits on m x and £cr p it does 
illustrate the qualitative effect of varying the properties 
of the halo model on the values of the WIMP parameters 
obtained from the data. A full likelihood analysis of the 
DAMA data has been carried out for models with bulk 
rotation in Ref. || , allowing us to check the reliability of 
Brhlik and Roszkowski's technique. 

Given the experimental difficulties of extracting a 
small annual variation from, possibly time dependent, 
backgrounds, concerns have been expressed about the in- 
terpretation of the earlier (1 and 2 year) DAMA data as 
evidence for a WIMP signal [^5| . Furthermore the Cryo- 
genic Dark Matter Search (CDMS) collaboration have 
recently released limits on the WIMP-Nucleon cross sec- 
tion which exclude, at 85% confidence, the entire DAMA 
3er allowed region |2(|. The DAMA collaboration have, 
however, performed a thorough analysis of the various 
sources of possible systematic errors [Effl and, since the 
DAMA and CDMS experiments use different target nu- 
clei, assumptions are required to perform a direct com- 
parison. In any case even if the interpretation of the 
DAMA annual modulation signal as WIMP scattering is 
eventually found to be erroneous, our results will still in- 
dicate the qualitative effect of these non-standard halo 
models on the analysis of a WIMP annual modulation 
signal. 

In Figs. ||- H we plot contours of k = 35, delineating 
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FIG. 5. The k — 35 contour for flattened halo models 
with q = 1.0,0.85,0.707 (solid, dashed and dotted lines re- 
spectively). 
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FIG. 6. The k — 35 contour for the standard halo model 
with a rot = 0.36, 0.5, 0.64 (dashed, solid and dotted lines re- 
spectively). The non-rotating halo corresponds to a ro t = 0.5. 



V. CONCLUSIONS 

In this paper we have examined the effect of an asym- 
metric halo velocity distribution, halo flattening and bulk 
halo rotation on the WIMP annual modulation signal. 
Whilst for each of the halo models the change in the 
mean signal is small, < 10% for the range of recoil en- 
ergies probed by the DAMA experiment, the amplitude 
of the annual modulation can change significantly. The 
magnitude of the change depends on the WIMP mass 
and recoil energy, as well the nature of the deviation of 
the halo model from the standard Maxwellian halo. With 
bulk rotation the change in the annual variation, relative 
to the Maxwellian halo model, is largest for small re- 
coil energies, whilst halo flattening produces the largest 
change for larger recoil energies. 

The range of WIMP masses consistent with the DAMA 
annual modulation signal is enlarged significantly, to 
roughly 30 < m x < 150 GeV at 3a, for each of the mod- 
els considered. The shape of the allowed region in the 
m x — £<T p plane is different for each model however. This 
indicates that the uncertainties in halo modelling have 
a significant effect on the WIMP mass determined from 
the annual modulation signal. Therefore more sophisti- 
cated halo models (such as those in Ref. |^|) need to be 
developed and used in the analysis of annual modulation 
data. 
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the region of m x -^a parameter space compatible with the 
DAMA annual modulation signal, for asymmetric, flat- 
tened and rotating halo models respectively. In Fig. |] 
we see that as the asymmetry of the velocity distribution 
is increased the allowed region is enlarged and moves to 
larger masses and slightly smaller interaction cross sec- 
tions. We saw in Fig. || that the change in the annual 
modulation signal for the flattened halo model is larger 
than that for the asymmetry velocity distribution. Con- 
sequentially the change in the allowed region for the flat- 
tened halo model is smaller, extending to larger masses 
and slightly larger interaction cross sections. In Fig. ^ 
we can see that counter-rotation (a rot < 0.5) contracts 
the allowed region and shifts it to smaller cross sections 
whilst co-rotation (et rot > 0.5) expands it and shifts it to 
larger cross sections. This is in good agreement with the 
results of Ref. [|| on the effect of bulk rotation, which 
were found via a full likelihood analysis of the DAMA 
data. 



APPENDIX A: INTERACTION CROSS SECTION 

The differential cross section for neutralino scattering 
off a target nucleus is dominated by the scalar and axial 
terms [§|l§Jl(J: 



da 



da' 



calar 



da 



axial 



dq 2 



dq 2 



(Al) 



where q = (niAm x ) / (mA+rn x )v is the momentum trans- 
fered. 

The scalar differential cross section, which arises due 
to Higgs boson and squark exchange, is given by {l9| 



da 



scalar 



i 



dq 2 



[Zf p + (A-Z)f n ] 2 F\q), (A2) 



where / p and / n are the effective neutralino couplings 
to the proton and the neutron respectively and F(q) is 
the scalar nuclear form factor. The form factor for Na 



5 



is usually taken to be equal to unity, whilst for I the 
Saxon- Woods form factor IpOtl 



qKi 



(A3) 



where i?i = y/ R 2 A - 5s 2 , R A = A 1/3 x 1.2fm and s = lfm, 
is used. 

The axial differential cross section, which arises due to 
Zq and squark exchange, is given by |hJ 



dcr 



axial 



dq 2 



-^A 2 J(J+l)S(g) : 



(A4) 



where J is the total angular momentum, S(q) is the spin 
form factor and A depends on the axial couplings of the 
neutralino to the quarks (see Ref. Jl9| for more details 
and explicit expressions). 

The differential event rate for a given detector can then 
be expressed as fl9|: 



cLR 
dE 



r 3/2 m 



T(E) 



(A5) 



x{[Zf p + (A- Z)f n f F 2 (q) + 8A 2 J( J + l)S(q)} . 
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